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Abstract 


Thu computation of accunte geoid undulations Is usually done combining 
potential coefficient Information and terrestrial gravity data in a cap sur- 
rounding the computation point. In doing this a spherical approximation is 
made that can cause errors that are investigated in this pnper. The equations 
dealing with ellipsoidal corrections developed by Lelgemann and by Moritz 
are used to develop a computational procedure considering the ellipsoid as a 
reference surface. Terms in the resulting expression for the geoid undulation 
are identified as ellipsoidal correction terms. These equations have been 
developed for the case where the Stokes function is used, and for the case 
where the modified Stdies funct ion is used. For a cap of 20° the correction 
can reach -33 cm. 

Ellipsoidal corrections were also computed for the Marsh/Chang geoids. 
These corrections reach -45 cm for a cap size of 20°. 

Global maps are given showing the distribution of the corrections so 
that more accurate geoid undulations can be found. 
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1. Ultroductlon 
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Cit'old undulations can be coin|)utcd from Stokes' equation Riven a 
Rlob:il set of Rravity anomalies. The S toko s^ equation is usually Riven with 
spherical approximation (Helskanen and Moritz, 1967, p. 94). The error of 
the spherical approximation has been investigated by several authors who have 
sought a solution of the boundary value problem with the ellipsoid as a refer- 
ence surface. Lelgemann ( 1970) developed techniques to compute corrections 
to the undulations computed from the Stokes' equation. He showed that the 
root mean square correction was i 0. 2 m with a maximum correction of about 
0.6m . A somewhat revised approach to this approximation problem is dis- 
cussed by Moritz (1980). 

Today, however, geoid undulations are not usually computed through 
a global Integration procedures. Instead potential coefficient information is 
combined with gravity awmaly Information in a cap surroudning the com- 
putation point (Happ and Hummel, 1975, Marsh and Chang, 1976, Rapp, 1980). 
For this procedure vre therefore cannot use directly the correction equations 
developed by lelgemann; instead we must develop new equations. 


\ 
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2. The Theory 


We start with the following representation of the earth's gravitational 
disturbing potential: 


(1) T(r, <p. 

X) = Z (-§-) t. (Cn, cosmx + Sn. sinmX) P„, (sin<p) 

Brra B = 0 

where: 


kM 

geocentric gravitational constant; 

r 

geocentric radius; 

a 

earth equatorial radius used in potential coefficient 
determ inations ; 

^n« * ®nn 

potential coefficients. (The Ca,o * C 4 ,o , Ce,o etc. are 
referenced to a defined reference ellipsoid.) 

Pn> (Sln^ , 

associated Legendre functions, as a function of geocentric 
latitude o . 


The disturbing potential can be evaluated on the ellipsoid (assume no external 
masses) by letting r = re the geocentric radius to a point on the ellipsoid. 

We now can compare this to the following representation of T on the ellipsoid 
(Moritz, 1980, p.320): 

® a 

(2) T(o,X) = 2, r (A„, cos mX + Bn, sinmA) Pb, (sln^ 

B=a i=.’0 

where An, , Bn, are coefficients to be determined; <f> is the geodetic latitude . 
We need to find the relationship between the A , B coefficients and the C , S . 
coefficients noting that at the same point, equation (1) and (2) must yield the 
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same potential. We write ( in a linear approximation). 

?^T 

(3) T((p,^, A, B, ) T(ip,X, A, B) + (p - ?5) 

<9(p 

= T(® ,X. A, B) + 


r r (An,cosmX+ Bn, slnmX) - 


a =a , = o 

Now from geometric geodesy we know; 

(4) (</J -0) -Ai e^ sin?>cos^ 

From Moritz (ibid, 39-46, ) we can write; 

(5) sln0 cos^o " * “ 2 .„ Pp-a, , + b„ P** Pmji,, 

d0 

where; 


(6) bn» 


-n (n - m + iwn - m + 2i 
(2n + l)(2n + 3) 

n^ - 3m^ + n 


'' ““ (2n + 3)(2n-l) 

c = in + IWn + m)(n + m - 1 ) 

(2n + l)(2n-l) 

We now use (4) and (5) in (3) to write; 

(7) T(^,X,A,B) = f f (An. cosmX+ Bn, sinmX) P„, (sin^J) 

D =a ■ =. 0 

+ e^f E (An.cosmX + Bn. sinmX) (an, Pn«„ + bn, Pb. + c„. P„_*,, ) 

n = a ,_o 

To simplify (7) we use the following relationships (ibid, eq. 39-48); 
f a., {g-} , f a„.,. {*■«■■} P.. 

S jAn,l. ® /AaJ8,a\ 

n=a n=* Bira,a 

Now we can write (7) as follows; 

«> n r 

(9) T((p,X, A, B) = E E L(An. + e®Kn.) cosmX 

B = 0 • = 0 _ 

+ (Bn, + e® Ln, ) sinmX J Pn, (sing) 


where; 


Kna A„_a,a + b,, A„, + ^ n+a, a A^^, 

Bna ~ ®'n-a,B Bn-3,, + bn, Bb. + Bn4a,, 
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(Note tlvU the summation on n has started fmm zero as a consequt'noe of 
the n'latlonships in {$).) Wc now equate (9) and (1) with r = q (litd 


<») ^ (^)" (s:;} 


f A„ . + e‘ 


I 


a. + e‘ 



The A„, , Bn, coefficients are a function of i|» since JTt is a function of 0 . 
in a spherical approximation (11) becomes 



We next consider the ellipsoidal form of the Stokes' equation as 
given by Moritz (ibid, eq. (49 -26)): 


(13) Ne 

Here we have assumed the solution of the geoid boundary value problem as 
opposed to the Molodensky problem as discussed by Moritz. In (13) we have; 
y, average value of normal gravity over the ellipsoid; N, approximate value 
of the geoid undulation. In addition (Moritz, Ibid. , 39-80): 


-A. 
4 ny 


// 


(Ag - e® Ag') S(;/))d(T+ e^ ( 


1 

T 


4 


sin^<(l^ 


N 


(14) Ag' = 



(Gn,cosmX+ Hn, slnmX) Pn,(sln^) 


where; 


(18) 

with 

(19) 






J 

■ 

~ X nB J ■*" ^»B In J 

Bn-a ,1 “n , 

•f 

t ft o J 

H n» ' 

= -3(n - 3)(n - m - l^n - ml 




2(2n-3)(2n-l) 



\ 

„ n'^ - 3m‘^n - 9n^ - 6m^ - lOn 

® + 

_a. 

A 

3(2n + 3)(2n - 1) 




_ -(3n + 5Wn + m + 2wn + m + 

JUl 


M n a 

2(2n + 5)(2n + 3) 




Although Moritz starts the summation from n = 2 in (14), the formality 
of the reduction leading to (14) indicates a starting summation from n = 0. 
We also have da= co8</5d(pdX (Lelgemann, 1970). 


Now consider the evaluation of (13) with the Ag values considered in a 
cap Oe surrounding the computation point, and in the remaining area of 
the sphere. We have 

, ifo ^ 2 TT ^TT 2 TT 

(20) Ne = 


n " TT .TT Z rr 

(Ag - e= Ag') S((^)) da + 4 ^ f J - e® Ag') S(0) da 


where 

(20a) A N = e® ( ^ - -|- sin® (p) N 


In (20) (^0 is the spherical cap radius surrounding the computation point. 
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We introduce the new function S (|/)) (Heiskanen and Moritz, 1967, p. 259): 

(o It 0 ^ il) < il>o 

(21) 

Equation (20) then becomes*. 

.^/.2n 


/•i(b/'2TT 

(22) Nt = T^Jo Jo 4^y J 


S ((/) ) d a + an 


We now write (Moritz, ibid, p.326): 
(23) Ag^Ag° + e^Ag» 
where 


wnvxv « ft 

(24) Ag° =• J] Y. (A,, cosmX+ B„sinmX) Pn_(sin<p) . 

«=a 1=0 " 

n -2' 


Using (23) for Lg^ in (22) we have: 
•00^ 2 n 


fiiior^v -n,2Ti 

(25) Nc = J Ag°S(i/))da+ an 

Foiicwing Heiskanen and Moritz (1967, p.259) we introduce the Molodecskil 
coefficients : 

(26) S(0) = I Q„P„(cob(/)) 

8=0 

SO that the second integral in equation (25) can be written as: 

/ TT ^2 TT ® 

J Ag° S^(<^))d(T= ■^^L^(2n + 1) Qn r r Ag° P.(cos^/)) da 

Jo Jo Jo 

- £ Q.igs 

2y n=a 

where Ag° would be obtained from (24). 

If we now use (24) for Ag° and use (11) for the coefficient relationships 
we can express (25) as: 

(28) Ne = 4 ^fJ <^K -e Ag’)S(0)da 

+ r Q (n - 1) £ (Ca,cosmX 

2rty t ' \ rj/ « = o 

+ Sa, sinmX) Pn,(8in^) 

a <*> n 

" oV I ■ 1) E <Kft-CosmX + La, sinmX)Pa, (sin<p) + AN 

The summation cu n is from 2 if we ignore the zero and first degree terms 
that have been eliminated through the Stdces' equation. Note that the correction 
terms involving K and L arise from the fact that the associated Legendre 


1 ... 


■ ■ ,* J-"**. A- 


functions arc evaluated at the geodetic latitude Instead of the geocentric lat- 
itude. If the evaluation Is done at the geocentric latitude the correction terms 
will not t>e necessary. We have* 

00 2 n 

(29) Ne f (Ag-e=’Ag’)S(0)da 

Jo Jo 

2 ^^ Q»(D-1) (~y i (C^.cosmX 
2 ryy n_a \ rt / ,_o 

+ SoiSinmX) Pn, (sln‘5)) + AN 
Now consider the Integral in (29) involving e^Ag' . We write: 


0<v,2n 


17 

^ n ^ cs 


n 2n 


,3 ^^-20 


Ag’ S(0 ) d a 


jT J Ag’ S(0) do - Ag’S(0)do 


The first integral on the right hand side of (30) can be written as (Moritz, ibid 
p, 426): 

(31) J Ag’S(0)do= f 

+ Hn, sinmX) Pn« (sln^) 


Note the summation starts from n = 2 because of properties of the Stokes' 
fimctlon. Recall equation (14), however, where Ag' contains zero and first 
degree terms. Ihtroducting the modified Stokes' function (equation 21), the 
second Integral on the right hand side of (30) can be written as: 

(32) ^ f" f ig'S(*)d(T= f Ag;<3. 

J^oJo ^ n=0 


Combining (31) and (32) we can write (30) as: 


0 if n <2 

((A)- 


Y (G„, cosmX + Ha, SinmX) P,, (sin^) 

-••rt 


We now can write our final result by re-writing (29) with (33) and (21): 

(34) Nt = J— / / igS(«da 

•o •' 0 

+ ■ 1) (^e)" t ®n. BinmX) Pa, (Sin5) 

B ® iL*~^ 

(‘^" " ^») I (Ga, cos mX + H„, SinmX) . P,, (sln(p) 
' 8=0 1=0 

+ e*(-^ - -|-sin^(0)N 
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The first two terms In (34) represent the computation of the undulation 
with cerUiln approxImathTns while the latter two terms represent the corrections 
needed to hilly ixifer the solution to an ellipsoidal reference surface. Thus we 
would design our eomfxitatlons in the following usage: 

(35) Nj " Ni + Na + A + AN 

where Ni -f Na represents the first two terms in (34) and ANx + AN re- 
present the ellipsoidal correction terms. In a later part of the paper we will 
discuss the numerical values of ANx + AN . 


We can look at two special cases of the above equations. First let 
4)q ~ which implies that no gravity anomalies are used in the computation. 
In this case ■- 2/(n - l)(n 2) so that (34) becomes: 

(36) Me = (^y ^ (C„cosmX+ S„,sinmX) P„(sin5) 

+ i-J - T ^ 


Kquation (36) is the same as: 


(37) Ne 


kM 
r z yz 


£» / \ B 

I L (Cni cosmX + S„, sinmX) Pn. (8in<P) 

n = a ^ ■ = Q 


since (Moritz, 1980, eq. (39-17)): 

(38) Yt = y(l - e® + ^ e® sln^<p) 


Equation (37) is the same as given in Ilapp (1967, eq. 7) for the computation 
of a "rigerous" gcoid undulation. We thus see the satisfactory reduction 
of the general case derived in this paper to the special case previously known. 


The second special case to consider is when tpo = 180® . Ihen Qn = 0 
so that (34) becomes: 

m Nt - r r igs«i)do 

•'0 ®'o 

” 1 “ 

~ L / n ~T\ T (G„»cosmX+ H„, sinmX) P„, (sine?) 

^ »=a( , = o 

+ e^ ( -J- - ^ s in N 

This result can be compared to the solution given by Lelganam (1970, eq.(3-6)) 
w'here numerical tests show that the two formulations yield the same cor- 
rections to about -1.6cm. We thus have an additional special case confirma- 
tion of our general formula. 
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3. Ellfpsotdal Corrections Using the Modified Stokes’ Equation 


TesUs described by Rapp (19B0) and Jekeli (1980) indicated a signi>- 
ficant improvement in geoid undulation determinations if the Stokes' func- 
tion is modified by subtracting t»So from S(^) . This procedure 

can be represented in our case by rc'Writing equation (13) in the following 
form: 


(40) Ne 


4 Try 
R 


f (Ag-e^Ag')(S(0)-So)da 

J J (^g - e® Ag') Soda 


+ ifyj J (^8 - e®Ag') S(0) da 

a-Oc 

+ AN 


Using a procedure followed before considering jekeli (1980, sections 2 and 3) 
we arrive at a result similar to equation (34)* Specifically: 

R f>^»2v 

(41) Ne - J Ag(S(i()) -Sc^da 

^ JiM. f (n - 1) (-f-) r (C„, cosmX + S„ sinmX) Pn, (sin5) 
2rcy„=a .to 

gS 08 _ a ^ 

+ rv" r (Ga,cosmX+ Hn, SinmX). Pn,(sin<P) 

^ n=Q »-0 

+ AN 

We have (Jekeli, 1980, eq. 65): 

(42) Qn(0o) = Qn (0o) + (Pn -1 (COS0p) - COS 00 P. (COS 0o)) ; n s 1 

EquaUon (41) is the same as (34) with two exceptions: S(0) - S© replaces S(0) 
and Q„ replaces Q.^ . Numerical tests of both equations will be described 
in the following section. 



4. Numerical Results 


We now will evaluate tie ellipsoidal correction terms AN^ and AN 
defined in equation (35) and (34) and the similar terms in e luation (41). Our 
starting potential coefficieits are those of GEMS (Lerch et al, 1979) taken to 
degree 20. The first step in the confutation is to find the A., , B.. coefficients 
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using equation (11^. 

We next deter mined ANi4 AN as given Inequation (35). This Is 
elllptlclty correction Uk bo added to a geold undulation computed from 
the first two terms on the right hand side of equation (34). This evaluation was 
done for 10®, 20® , and 180° and the results are shown In Figures 1, 2, 
and 3. The maximum correction and root mean square correction for each of 
of these oases Is • (-26cm, dOcm, 0 =10®); (-33cm, ±10 cm, =20°); (-59cm, 

±18 cm, ;^)~180^). 

Simlliar computation were carried ovt when using the modified Stokes' 
equation. These results are shown for 0= 10® in Figure 4 and for tl)= 20° In 
Figure 5. The maximum and the RMS correction for 10° Is (-21 cm, ±5 cm) 
and for 0 = 20° It is (-27 cm, ±6 cm). Examination of Ihe corresponding figures 
indicate that the correction for the modified Stokes' integral are somewhat smaller 
overall than the case with the regular Stokes' function. 

The corrections are generally small and below the current acccxracy of 
the data w ith caps of 10® or 20®. However as more precise computations 
are carried out in the future, these correction? should be taken into account. 


5. The Zero and First Degree Ccrrection Problem 


In carrying out the derivation of several of the previous equations, sum- 
mations were started from 2 instead of 0 by convention or because the Stokes' 
equation removes '/sero and first degree terms in a global Integration. However 
the use of the relationships in equation (8) does Introduce zero and first degree 
terms that need to be considered. This problem has been discussed by belgemann 
(1970) who assumed the following form of the disturbing potential: 


(43) T(r,9,X) = I Aa,oRB,o(5, X) + Aa^ Ra^ (0,X) ± Ba ,9 Sa ,8 (S|X) ) 

^3fi Rs,o(5»X) + Aa,i Rg,i (5 ,X) + Aa ,3 Ra,a (5,X) 

+ Ag^a Ra,a ( 0 , X) + Ba,i Ro,i (5 , X) ± Bg^a S 3 , a (S . X) 
± Ba,3 S 3,3 (5, X) 1 + P 4 ,1 (5, X) ] 

where; 9^ = 90° - 0 

R„, = cos mXp^, 

S„, = slnmXPn* 


|A„| 

• - kMa" (^"•1 


nm 


Note that tlie are not the same as given in equation (2). Equation 

(43) represents a reasonable low degree model of the disturbing potential but 
it is not meant to be a complete model. 
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Undulation Correction When A Global Set of Anomalies 
are Used, ^ = 180° ; (Contour Interval = 5 cm) 












Lolgemann (1970) shows that the error Introduced by the neglect of the 
zero and first degree terms when developing the ellipticity correction for the 
global Stokes' integration ((^ = 180®) Is; 


(44) ANoa 


c* 

y 



12 


‘ a, 0 


Ri, 




Ru + ^3,1 Sii ) ) 


where; 


— k M 

Ibn.J a \s„. 


Using the coefficients of GEM9 we have evaluated equation (44) with the results 
shown in Figure 6. As is obvious from (44) this correction is very long 
wavelength. The magnitude is quite small with the maximum correction being 
l7cm. We would expect that the correction for the small cap sizes used in 
practice would be considerably smaller than this as was seen for the usual 
ellipticity correction. Therefore we will not pursue the derivation of this 
correction term for the cap case. 


6* Ellipticity Corrections for the Marsh-Chang Geoid 

For the past several years Marsh and Chang (1976, 1978) have computed 
detailed geoid undulations combining potential coefficient Information and ter- 
restrial gravity data. The method used by them is called Method A in Rapp and 
Rummel (1975) or Method 1 in Rummel and Rapp (1976). The specific equations 
used by Marsh and Chang in their recent papers are as follows: 


(45) N(9,X) = N. +N« 


where; ^ „ 

(46) Nj = R £ (CBaCosmX+ Sn,sinmX) (sin^ 
n=a 9 ~o 


(47) Na = ^ j I (Ag-Ag.) S{ii))do 


JJc 

n n 


(48) Ag, = y r ~ sinmX) Pn, (sin<p) 


n~3 , = o 


Here ii is the maximum degree used with the potential coefficients. In 
practice the Integration in the cap has taken place us ing just 1° x 1° anom- 
alies or these anomalies In conjunction with smaller block sizes such as 5' x 5', 
and 15' x 15' . The integration cap has been 10°, 20° or presumably o° 
when insufficient gravity data is present. In view of our previous discussions 
we now are interested in the geoid error caused by the spherical approximation 
in equation (46) and (47). 


For convenience we introduce A^ as follows; 
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(44) Ab = |[ (C„cosmX+ Sn, slnmX) P„ (sln<3) 

■ = 0 

The ellipsoidal error in the Marsh/Chang geold will be; 


(50) AN 44 /C = Nc (eq. 34) - N(eq. 45) 

In order to reduce (50) the following equality cau be used 

(51) 4 ^ / / ^ “ -f A. 


An erroneous form of (51) was given in Rapp and Rummel (1975, eq. 33). 
Using (34), (45), and (51), equation (50) can be written as: 


(52) “ 57^ E Q.(»-l) (-^)' A. 

- A, 


.3 « 


+ T. 


X„) f (GB,cosmX 


n = 0 B aO 

+ Hbb sinmX) Pb, (sin^ + AN 


We have evaluated (52) with the GEM9 potential coefficients (fi = 20) for 
|/)= 0°, 10°, and 2(f . These results are shown in Figure 7, 8, and 9. The 
maximum correction and root mean square correction for each of these cases 
is (101cm, ±27 cm, (^>=0°), (44cm, ±16cm, (/)=10°), (-45cm, tl4cm, 

l/)= 20° ). 


7. Summary and Conclusions 

This paper has developed the formulas needed to compute the correction 
for geo id undulation compiiations made from the combination of potential coef- 
ficient information and terrestrial gravity data. The first procedure developed 
the formulas needed for the precise computation of the geold considering the 
ellipsoid as a reference surface and using the usual Stdees' equation. The 
corrections are a function of the cap within which gravity data is used. For 
a cap size of 20° the maximum correction was -35cm. 

Another case was considered with the use of the modified Stokes' function 
In this case the maximum correction for (/)= 20° was -27cm. 

The third case considered was that for the elliptic ity corrections for Ihe 
Marsh/Chang geoid. If a cap of 0° was used the maximum correction was 
101 cm; if the cap was 20°, the maximum correction was -45cm. 

All the corrections have been computed using the GEM9 potential coef- 
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ficients taken to- degree 20. Errors in Uicse coefficients or the use of ad- 
ditional higher degree terms should not significantly effect these results. 

As can be seen from Uie various maps, these correction tt'rms are 
fairly long wavelength. Therefore in some applications worldng with alti- 
meter data and gravimetric gcoids, the correction could appear as a constant 
difference. In some cases, for example, in examining the difference between 
the sea surface and the gravimetric geoid of Marsh/Chang a net correction, 
across the Pacific Ocean, of about 35cm should be made. If, in the future, 
we are to determine highly accurate geoids from potential coefficients and 
terrestrial gravity data, the corrections or problem formulation given In 
this paper should be used. 
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